library(sf)
ehime_map0 <- sf::read_sf("N03-20220101_38_GML/N03-22_38_220101.shp", options = "ENCODING=CP932")
office_map0 <- sf::read_sf("P34-14_38_GML/P34-14_38.shp", options = "ENCODING=CP932")単純な巡回セールスマン問題
この春,たかしくんは愛媛県庁に入庁した。 これから市町村役場の方々と一緒に仕事をする機会も増えるだろう。 さっそくみんなに顔を覚えてもらうために,愛媛県内にある市区町村役場等の施設にあいさつ回りに行こうと考えた。 それでは,どのようなルートで訪問すればよいだろうか。
1 データの取得
まず,愛媛県の地図データと市町村役場等の施設の位置情報に関するデータをダウンロードすることから始めよう。 国土交通省国土数値情報ダウンロードサイトから,行政区域データと市区町村役場データをダウンロードする。 ファイル名はそれぞれ,N03-20220101_38_GML.zipとP34-14_38_GML.zipである。
ダウンロードしたZIPファイルを解凍し,Rでシェープファイルを読み込む。
このようにしてもよいが,kokudosuuchiを使うとZIPファイルを解凍する必要がなく,便利である。
library(kokudosuuchi)
ehime_map <- kokudosuuchi::readKSJData("N03-20220101_38_GML.zip")
office_map <- kokudosuuchi::readKSJData("P34-14_38_GML.zip") |>
translateKSJData()それぞれの方法で読み込まれるデータはまったく同じではないことに注意が必要である。 後者はリストのリストとして読み込まれる。
is.list(ehime_map0[[1]])[1] FALSE
is.list(ehime_map[[1]])[1] TRUE
市町村名をfactorにしておくと,図を描いたときの凡例の順序がおかしくならずにすむ。
ehime_map[[1]]$N03_004 <- factor(ehime_map[[1]]$N03_004, levels = unique(ehime_map[[1]]$N03_004))1.1 データの修正
先ほどダウンロードした愛媛の令和4年の行政区域データには,市区町村名(N03_004の列)がNAとなっているレコードが存在する。
table(is.na(ehime_map[[1]]$N03_004))
FALSE TRUE
4906 1
この部分を残したままだと,以下の図の凡例にNAという名前の市町村が現れる。 これをどう処理するかは考える必要があるが,たかしくんは強引にもレコードごと削除した。
ehime_map[[1]] <- ehime_map[[1]][-which(is.na(ehime_map[[1]]$N03_004)), ]2 地図のプロット
Rにデータを読み込めたので,さっそく地図をプロットしてみる。
# plot(ehime_map[[1]])# 処理に非常に時間がかかるため,避けた方がよい
# plot(sf::st_geometry(ehime_map[[1]]))
library(rmapshaper)
ehime_map_shrink = rmapshaper::ms_simplify(ehime_map[[1]], keep = 0.001, keep_shapes = TRUE)
par(mar=c(0, 0, 0, 0))
plot(sf::st_geometry(ehime_map_shrink))
plot(sf::st_geometry(office_map[[1]]), pch = 16, col = "orange", add = TRUE)
ぜんぜん悪くない。 ただし,たかしくんはggplot2を用いた図の方が好みである。
library(ggplot2)
# カラーパレット
library(viridis)
ggplot(ehime_map[[1]]) +
geom_sf(aes(fill = N03_004)) +
geom_sf(data = office_map[[1]], colour = "orange") +
scale_fill_viridis_d(option = "mako", direction = 1) +
guides(fill = guide_legend(title = "市町村名")) +
theme_void() +
theme(text = element_text(family = "HiraKakuProN-W3"))
たかしくんはここで重大な問題に気がついた。 この地図は動かせないし,拡大もできない。 Google マップに慣れ親しんでいるたかしくんにとって,地図は動かせた方が便利である。 そこで,Leafletを使うことにする。
library(leaflet)
micanfIcon <- makeIcon(
iconUrl = "https://www.pref.ehime.jp/h12200/mican-kanzume/images/img_dance02.png",
# iconUrl = "https://www.pref.ehime.jp/h12200/mican-kanzume/images/img_dance03.png",
iconWidth = 40, iconHeight = 40,
iconAnchorX = 20, iconAnchorY = 20,
)
leaflet::leaflet(office_map[[1]]) |>
addTiles() |>
addMarkers(icon = micanfIcon)みきゃんだ。
library(mapview)
mapview::mapview(ehime_map[[1]], zcol = "N03_004", layer.name = "市町村名") +
mapview(office_map[[1]], col.regions = "orange", legend = FALSE)